function lambda_hat=find_lambda_hat(w,w_grid,h,h_density,tau,kernel_type,order)

N_grid=size(w_grid,1);

	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
	% LIMITS
	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
delta=(w_grid-w_grid')/h;
x_ub=min(1,1-delta);
x_lb=max(-1,-1-delta);

check=(x_ub<x_lb);
x_lb(check)=x_ub(check);

Cmat=NaN*ones(N_grid,N_grid);

	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
	% C MATRIX
	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
for ii1=1:N_grid
	fhat=density_X(w,w_grid(ii1),kernel_type,h_density);
	for ii2=1:N_grid
		Cmat(ii1,ii2)=tau*(1-tau)*fhat* (fct_C_inner(x_ub(ii1,ii2),delta(ii1,ii2),kernel_type) - fct_C_inner(x_lb(ii1,ii2),delta(ii1,ii2),kernel_type));
	end
end

	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
	% COMPUTE EIGENVALUES
	%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
[~,lambda_mat] = eig(Cmat);

lambda_hat=sort((1/N_grid)*diag(lambda_mat),'descend');
lambda_hat=lambda_hat(1:order);

% 





